function positions = Lattice(Nx,Ny)
%LATTICE collect coordinates of each atom
N1=Nx*Ny; % on-lattice
N2=(Nx-1)*(Ny-1); % off-lattice
N=N1+N2;

positions=cell(1,N);
xi=1;
yi=sqrt(3)/2;

% set coordinates of each point
for ii=1:Ny
    if mod(ii,2)==1 % odd
        for jj=1:Nx
            positions{(ii-1)*Nx+jj}=[(jj-1)*xi -(ii-1)*yi];
        end
    else
        for jj=1:Nx
            positions{(ii-1)*Nx+jj}=[(jj-1+1/2)*xi -(ii-1)*yi];
        end
    end
end
%
for ii=1:Ny-1 
    if mod(ii,2)==1 % odd
        for jj=1:Nx-1
            positions{N1+(ii-1)*(Nx-1)+jj}=[(jj-1+1/2)*xi -(ii-2/3)*yi];
        end
    else
        for jj=1:Nx-1
            positions{N1+(ii-1)*(Nx-1)+jj}=[(jj-1+1)*xi -(ii-2/3)*yi];
        end
    end
end


end

